Chapter 1: Drude Model¶
Drude model considers electrons small solid particles traveling through lattices as a stationary array of heavier and bigger atoms. The model use "mean free time" as a measure of the frequency an electron gets scattered. This jupyter notebook makes use of this concept and uses stochastic method to simulate electrons traveling through a 2D material of finite width but infinite length.¶
In [1]:
import matplotlib.pyplot as plt
import numpy as np
import ipywidgets as ipw
from numpy import random
from ipywidgets import interact
from matplotlib import animation
from IPython.display import HTML
%matplotlib widget
First we setup the simulation parameters such as the external electric field, number of electrons (default to 1) as well as the pseudo-time steps to take. Also, we define the various physical constants such as elementary electron charge, (bare) electron mass, and mean free time ($\tau$). For the simplicity of simulation, we also define a constant called $V_0$ as an average restarting velocity after a collison.¶
In [2]:
e = -1.602E-19; # Electron charge (C)
m = 9.109E-31; # (Bare) electron mass (kg)
v0 = 3; # Average restarting velocity (m/s)
steps = 1000; # number of updates to simulate
Now we setup up four cases for easy switching in between with the first case corresponding to:¶
1. Only Electric field with long mean free time.
2. Only Magnetic field with long mean free time.
3. Having both electric and magnetic fields with short mean free time.
4. Duplicate of case 3 as a playground allowing adjustment to the variables.In [3]:
# Set time step
global Ey,delta_T
delta_T = 1E-12; # time step in seconds
case = 3;
if case == 1:
num = 1; # number of electrons to simulate
tau = 2000; # Mean free time (s)
Ex = -0.2; # X component of the external electric field (V/m)
Ey = 0; # Y component of the external electric field (V/m)
Bz = 0; # Z component of the external magnetic field (T)
yBnd = 3*num*5E-11; # physical boundry of the conductor in y direction
elif case == 2:
num = 1; # number of electrons to simulate
tau = 8; # Mean free time (ps)
Ex = 0; # X component of the external electric field (V/m)
Ey = 0; # Y component of the external electric field (V/m)
Bz = 0.2; # Z component of the external magnetic field (T)
yBnd = 15*num*5E-11; # physical boundry of the conductor in y direction
elif case == 3:
num = 1; # number of electrons to simulate
tau = 20; # Mean free time (ps)
Ex = -1; # X component of the external electric field (V/m)
Ey = 0; # Y component of the external electric field (V/m)
Bz = 0.3; # Z component of the external magnetic field (T)
yBnd = 15*num*5E-11; # physical boundry of the conductor in y direction
elif case == 4:
num = 6; # number of electrons to simulate
tau = 20; # Mean free time (ps)
Ex = -1; # X component of the external electric field (V/m)
Ey = 0; # Y component of the external electric field (V/m)
Bz = 0.3; # Z component of the external magnetic field (T)
yBnd = 15*num*5E-11; # physical boundry of the conductor in y direction
else: # Default setting
num = 1; # number of electrons to simulate
tau = 10; # Mean free time (ps)
Ex = -0.2; # X component of the external electric field (V/m)
Ey = 0; # Y component of the external electric field (V/m)
Bz = 1; # Z component of the external magnetic field (T)
yBnd = 3*num*5E-11; # physical boundry of the conductor in y direction
Next, base on the choice of the scenario, we initiate some containers to store variables in the simulation¶
In [4]:
global x, y, vxt, vyt, ayt
x = np.zeros(num);
x = np.expand_dims(x,axis=0);
#y = random.rand(num)*yBnd-yBnd/2;
y = np.zeros(num);
y = np.expand_dims(y,axis=0);
vx = np.ones(num); # Initial velocity in x direction
vy = np.ones(num); # Initial velocity in y direction
vxt = vx; # Instantaneous velocity in x-direction
vxt = np.expand_dims(vxt,axis=0);
vyt = vy; # Instantaneous velocity in y-direction
vyt = np.expand_dims(vyt,axis=0);
ayt = np.zeros(num); # Acceleration in y-direction
ayt = np.expand_dims(ayt,axis=0);
colors = []; # Color assignment to trace different electrons
for i in range(num):
colors.append('#%06X' % random.randint(0, 0xFFFFFF))
Next we create the figure objects, and the function to call between each frame for simulation of each time step. In each step, we compute and plot the average velocity in x direction over the steps simulated so far on top of the instantaneous velocities. This mean velocity in x direction eventually should converges to a called drift velocity, which is an indirect measure of the conductivity/resistivity of materiasl in Drude model¶
In [5]:
# Create figure objects
figs = plt.figure();
ax1 = plt.subplot(1,9,(1,3));
plt.rcParams['figure.figsize'] = [10, 4];
init, = ax1.plot([],[],'ok', markersize = 5, zorder = 10);
end = ax1.scatter(x[0,:],y[0,:], c = colors[:], s = 30, marker = 'o', zorder = 10);
ax1.axhline(y = -yBnd, color = 'k');
ax1.axhline(y = yBnd, color = 'k');
ax1.set_ylim(-1.2*yBnd,1.2*yBnd);
ax1.set_xlabel('x-Distance (m)');
ax1.set_ylabel('y-Distance (m)');
ax2 = plt.subplot(2,9,(6,9));
ax2.set_xlabel('Time (s)');
ax2.set_ylabel('Acceleration (m/s^2)');
ax2.set_xlim(0,steps*delta_T);
ax3 = plt.subplot(2,9,(15,18))
ax3.set_xlabel('Time (s)');
ax3.set_ylabel('Velocity (m/s)');
ax3.set_xlim(0,steps*delta_T);
plt.show()
# Set the electrical field in y direction
def set_field(E):
global Ey
if(E==1):
Ey = -Bz*Ex*e*tau/m/10**12; # Y component of the external electric field (V/m)
else:
Ey = 0;
return
# Initialize the frame
def initf():
init.set_data(x[0,:],y[0,:]);
return init,
# Animation function which updates figure data. This is called sequentially
def update(step,Ex,vx,vy,tau,v0):
global Ey,x,y,vxt,vyt,ayt,delta_T
newX = np.zeros(num);
newY = np.zeros(num);
newVt = np.zeros(num);
ay = np.zeros(num);
for nn in range(num):
vx_old = vx[nn]; # Use a temporary variable to store the previous value
vy_old = vy[nn]; # Use a temporary variable to store the previous value
vx[nn] = vx[nn] + e*(Ex - Bz*vy_old)*delta_T/m; # Update the velocity
vy[nn] = vy[nn] + e*(Ey + Bz*vx_old)*delta_T/m; # Update the velocity
ay[nn] = (Ey + Bz*vx_old)*e/m; # acceleration in y-direction
newX[nn] = x[step-1][nn] + vx[nn]*delta_T; # Update the coordinate
newY[nn] = y[step-1][nn] + vy[nn]*delta_T; # Update the coordinate
if(random.rand() < 1/tau):
theta = random.rand()*2*np.pi; # pick a random angle in 2D
vec = random.rand()*v0;
vx[nn] = np.cos(theta)*vec;
vy[nn] = np.sin(theta)*vec;
# Restrict movements of the electrons in y-direction by a hard and elastic physical boundary
if(newY[nn] > yBnd):
newY[nn] = yBnd;
vy[nn] = -np.absolute(vy[nn]);
elif(newY[nn] < -yBnd):
newY[nn] = -yBnd;
vy[nn] = np.absolute(vy[nn]);
# newVt[nn] = np.linalg.norm([vx[nn],vy[nn]]);
x = np.concatenate((x, newX[None,:]), axis=0); # Append the new coordinate
y = np.concatenate((y, newY[None,:]), axis=0); # Append the new coordinate
vxt = np.vstack([vxt,vx]); # Append the average of x-velocity of the current step
vyt = np.vstack([vyt,vy]); # Append the average of y-velocity of the current step
ayt = np.vstack([ayt,ay]); # Append the average of y-acceleration of the current step
for ii in range(num):
# ax1.plot(x[len(x)-1,ii], y[len(y)-1,ii], '.', markersize = 3, color = colors[ii], zorder = 0);
# ax3.plot([xx*delta_T for xx in range(step+1)], vxt[:,ii], color = colors[ii], zorder = 0);
ax1.plot(x[len(x)-1,ii], y[len(y)-1,ii], '.', markersize = 3, color = colors[ii], zorder = 0);
ax2.plot([xx*delta_T for xx in range(step+1)], ayt[:,ii], color = colors[ii], zorder = 0);
ax3.plot([xx*delta_T for xx in range(step+1)], vxt[:,ii], color = colors[ii], zorder = 0);
end.set_offsets(np.vstack((x[len(x)-1,:],y[len(y)-1,:])).T);
ax2.plot(step*delta_T,np.mean(ayt),'.',markersize = 3,color = 'k');
lgda =np.append(['Acc_y']*num,'mean(Acc_y)');
ax2.legend(lgda,loc='upper right');
ax3.plot(step*delta_T,np.mean(vxt),'.',markersize = 6,color = 'k', zorder = 10);
lgd = np.append(['Vx']*num,'mean(Vx)');
ax3.legend(lgd,loc='upper right');
return end,
# Call the animator. blit=True means only re-draw the parts that have changed.
anim = animation.FuncAnimation(figs, update, init_func=initf, fargs = (Ex,vx,vy,tau,v0), frames=range(1,steps+1), interval=100, blit = True, repeat = False)
interact(set_field, E=ipw.FloatSlider(min=0, max=1, step=1, value=0, description = 'Ey (off/on)'));
interactive(children=(FloatSlider(value=0.0, description='Ey (off/on)', max=1.0, step=1.0), Output()), _dom_cl…
Now as an exercise, you may try to find values in real materials for "Mean Free Time" and appropriate time step to run simulate again and see if you could find the right "Average Velocity" for electrons traveling in these materials.¶
Thank you for using the notebook, please provide us with some feedbacks at your convenience by emailing to (yikai.yang@epfl.ch), (oleg.malanyuk@epfl.ch)¶
In [ ]: